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FINITE-DIFFERENCE  SIMULATIONS  OF  RAYLEIGH  WAVE 
SCATTERING  BY  2-D  ROUGH  TOPOGRAPHY 

SUMMARY 

Rayleigh  waves  normally  incident  upon  2-D  simple  or  rough  topographic  structures  are 
simulated  by  the  linear  finite-difference  method  to  study  the  attenuation,  transmission,  and 
reflection  of  Rayleigh  waves  and  to  measure  the  Rayleigh-to-P  and  -SV  bodywave  conversion. 
For  simple  ramp  structures,  transmission,  reflection,  and  scattering  depend  on  the  sign  of 
change  of  slope  of  the  topographic  feature,  as  well  as  the  ratio  of  the  ramp  height  to  the 
wavelength,  h/X.  Simple  ramp  structures  produce  back-scattered  bodywaves  for  h  >X,  and 
forward-scattered  bodywaves  for  h<\.  The  radiation  patterns  of  P  and  S  bodywaves  are 
roughly  consistent  with  the  model  of  equivalent  point  forces  along  the  free  surface.  More 
complicated  topographic  features  generated  by  random  Markov  sequences  have  been  character¬ 
ized  by  the  Rayleigh-wave  spatial  Q(f)-  As  expected,  rougher  topography  attenuates  Rayleigh 
waves  more  than  smooth  topography.  P  and  S  amplitudes  ratios  are  consistent  with  radiation 
from  equivalent  point  forces  near  the  surface,  but  the  distribution  of  slownesses  generated  is 
greater  than  from  the  simple  ramp  structures.  Reflection  of  Rayleigh  waves  by  topographic 
slopes  and  by  random  topography  is  an  inefficient  process  and  the  bulk  of  the  energy  that  is 
not  transmitted  as  Rayleigh  waves  is  converted  to  bodywaves.  Fundamental  Rayleigh-to-Lg 
scattering  and  generation  of  teleseismic  P  coda  by  short  period  Rayleigh  should  be  observable. 
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INTRODUCTION 

With  some  simple  simulations  we  explore  the  stability  of  planar  Rayleigh  waves  incident 
upon  simple  inclined  surfaces  and  random  rough  topographic  relief.  We  examine  the  reflected 
and  transmitted  Rayleigh  waves  as  well  as  the  bodywave  field  at  depth  as  a  measure  of  the 
simulated  far-field  P  and  S  coda.  In  this  way  we  characterize  rough  topography  by  a  spatial  Q 
for  the  transmitted  Rayleigh  waves  and  examine  the  partitioning  of  energy  into  bodywaves. 

Since  most  previous  research  concerned  with  the  scattering  of  Rayleigh  waves  has  dealt 
with  the  simplest  geometric  configurations  (quarter-space,  three-quarter  space,  a  trench  or  a 
vertical  step)  where  the  boundaries  are  always  parallel  or  perpendicular  to  the  grid  axis,  we 
have  concentrated  on  the  inclined  surface  and  more  general  topographic  profiles.  Topographic 
scattering  through  a  right  angle  step  was  addressed  by  deBremaecker  (1958)  and  continues  as  a 
topic  of  interest.  Numerous  workers  have  examined  the  incident  Rayleigh  or  Love  wave  upon 
welded  quarter  spaces  (see  for  example  Mai  and  Knopoff,  1965;  McGarr  and  Alsop,  1967; 
Drake,  1972;  Munasinghe  and  Famell,  1973;  Tuan  and  Li,  1974;  Martel  et  al,  1977;  Sabina 
and  Willis,  1977;  Chen  and  Alsop,  1979;  Fuyuki  and  Matsumoto,  1980;  Drake  and  Bolt,  1980; 
Momoi,  1980,  1981,  1982;  Fuyuki  and  Nakano,  1984;  and  Levander,  1985).  Recent  work  has 
moved  to  more  complicated  geometries,  3-D  as  opposed  to  2-D,  or  to  more  complicated  struc¬ 
tures.  Ultrasonic  modeling  was  used  by  Bullit  and  Toksoz  (1985)  to  examine  obliquely 
incident  Rayleigh  wave  transmission  and  reflection  by  sloping  topography.  An  important  point 
made  by  Bullit  and  Toksoz  was  that  the  transmission  coefficient  is  not  phaseless  and  that  there 
are  phase  shifts  induced  by  the  scattering.  In  a  linearized  theory  of  topographic  scattering  of 
surface  waves,  Sneider  (1986)  points  out  that  the  scattering  loss  is  intimately  related  to  the 
phase  velocity  perturbation  induced  by  the  topographic  structure. 
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The  problem  of  Rayleigh  waves  scattered  by  topography  may  have  several  applications. 
Short-period  Rayleigh  waves  are  sometimes  observed  from  shallow  earthquakes  as  well  as 
from  underground  explosions.  The  scattering  of  such  short-period  Rayleigh  waves  by  topogra¬ 
phy  has  been  proposed  as  a  mechanism  for  the  generation  of  teleseismic  P  coda  by  Greenfield 
(1971).  Gupta  et  at.  (1986)  have  proposed  that  such  P  coda  may  be  diagnostic  of  the  near 
source  medium.  P  coda  observed  at  arrays  is  often  partly  composed  of  Rayleigh  waves  gen¬ 
erated  by  P  waves  steeply  incident  upon  shallow  heterogeneity  or  topography.  Key  (1967, 
1968)  identified  such  Rayleigh  waves  within  the  P  coda  at  the  EKA  array  (Eskdalemuir,  Scot¬ 
land)  for  certain  azimuths  of  arrival.  Dainty  (1986)  shows  by  use  of  frequency-wavenumber 
analysis  that  P  coda  at  the  small  NORESS  array  (southern  Norway)  is  often  composed  of  fun¬ 
damental  short-period  Rayleigh  waves.  Both  the  coherency  of  such  coda  and  its  persistence 
are  closely  related  to  the  scattering  of  the  Rayleigh  waves  by  topography  and  shallow  hetero¬ 
geneity.  Aki  and  Chouet  (1 975)  proposed  three  simple  models  for  the  generation  and  propaga¬ 
tion  of  short-period  coda  from  local  earthquakes.  They  considered  the  case  of  scattered  surface 
waves  as  well  as  scattered  bodywaves.  Although  the  coda  of  local  earthquakes  appears  to  be 
dominated  by  S  waves,  it  is  clear  that  trapped  modes  and  surface  waves  are  also  present  (Phil¬ 
lips  and  Aki,  1986).  A  portion  of  the  regional  Lg  phase  and  its  coda  is  composed  of  funda¬ 
mental  and  higher  mode  Rayleigh  waves  trapped  in  the  crustal  waveguide.  Der  eta l  (1984) 
propose  that  a  significant  portion  of  the  Lg  coda  is  composed  of  short-period  fundamental  Ray¬ 
leigh  waves  that  concentrate  energy  near  the  surface  and  are  multiply  scattered,  contributing  to 
the  long  Lg  coda  observed  on  regional  seismograms.  Baumgardt  (1985)  proposes  that  Lg  from 
Eastern  Kazakh  incident  upon  the  Urals  produces  observable  teleseismic  P  coda  at  NORSAR. 
All  of  these  problems  have  in  common  a  need  to  estimate  the  effects  of  scattering  on  short- 
period  surface  waves  due  to  shallow  heterogeneity  and  topography. 
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In  this  paper  we  address  the  problem  of  transmission,  reflection,  and  conversion  of  funda¬ 
mental  Rayleigh  waves  by  topography  on  a  homogeneous  halfspace  using  two-dimensional 
linear  elastodynamic  finite  differences.  We  first  examine  the  problem  of  transmission  through 
simple  slopes  and  the  partitioning  of  energy  by  the  scattering  process.  It  is  demonstrated  as 
expected  that  larger  scale  variations  in  topography  produce  larger  coefficients  of  reflection  and 
smaller  transmission  coefficients.  Reflection  is  always  observed  to  be  an  inefficient  process. 
We  have  also  observed  that  the  reflection,  transmission  and  diffraction  are  directional  depen¬ 
dent.  Both  directionality  and  frequency  dependence  can  be  seen  in  the  converted  P  and  S 
bodywaves  scattered  from  the  topography.  Next  we  examine  the  effective  spatial  Q  produced 
by  random  topography.  Attenuation  of  the  short-period  Rayleigh  wave  is  found  to  be  strong 
by  moderate  rough  topography,  and  Rayleigh-to-bodywave  conversion  is  seen  to  dominate  the 
attenuation  process. 
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THE  FINITE-DIFFERENCE  SCHEME 

The  2-D  linear  finite-difference  (LFD)  scheme  popularized  by  Kelly  et  al.  (1976)  is  used 
with  the  Clayton  and  Engquist  (1977)  (Emerman  and  Stephen,  1983)  absorbing  boundary  con¬ 
ditions  on  both  the  sides  and  the  bottom  of  the  grid.  Even  with  the  absorbing  boundary  condi¬ 
tions,  the  simulations  are  run  with  large  enough  grids  such  that  for  the  duration  of  the  simula¬ 
tions  reflections  from  the  sides  do  not  seriously  contaminate  the  results.  The  ffee-surface 
boundary  conditions  follow  from  the  work  of  Jih  et  al  (1986)  which  allows  for  the  introduc¬ 
tion  of  arbitrary  polygonal  topography  into  the  wave  propagation  simulations  and  has  been 
applied  by  McLaughlin  and  Jih  (1986)  for  the  modeling  of  teleseismic  P  waves  from  explosive 
line  sources  under  realistic  topographic  surfaces.  These  boundary  conditions  are  extensions  to 
those  of  Ilan  (1975  and  1977).  Recently  a  number  of  techniques  have  been  pursued  in  an 
effort  to  improve  the  computational  performance  of  LFD  solutions  to  wave  equations.  These 
include  variable  grids,  higher  order  schemes  or  implicit  rather  than  explicit  methods.  We  have 
found  that  none  of  these  methods  provide  significant  improvement  over  the  traditional  explicit 
second  order  centered  LFD  schemes  described  by  Alterman,  Boore  and  Kelly  when  rough 
topographic  profiles  are  used  at  the  free  surface  (Jih  et  al.  1986). 

In  the  present  experiments,  the  medium  has  P-  and  S-wave  velocities  of  5.0  and  2.96 
km/scc  respectively,  corresponding  to  a  Poisson’s  ratio  of  0.35  and  a  Rayleigh  wave  velocity 
of  2.71  km/sec.  The  formulation  generating  the  incident  pure  Rayleigh  wave  packet  is  analo¬ 
gous  to  Boore  (1970)  or  Munasinghe  and  Famell  (1973)  which  has  been  used  extensively  by 
Martel  et  al.  (1977),  Fuyuki  and  Matsumoto  (1980),  Fuyuki  and  Nakano  (1984),  and  Levander 
(1985).  The  initial  waveform  is  the  fundamental  mode  Rayleigh  wave  in  the  homogeneous 
half  space  with  a  Ricker  wavelet  shape  (Figure  1).  This  method  is  superior  to  another 
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approach  which  shakes  several  mesh  points  on  the  free  surface  or  the  grid  edge  to  generate  the 
incident  Rayleigh  wave  as  well  as  undesirable  bodywaves  (Wojcik  and  Mak,  1985).  To  reduce 
the  grid  dispersion,  we  have  oversampled  the  incident  wave  packet  by  choosing  the  dominant 


wavelength  to  be  3.2  km  =  32  Ax,  where  Ax  =0.1  km  is  the  grid  spacing. 
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NUMERICAL  EXPERIMENTS:  SIMPLE  INCLINED  SURFACES 


Experiments  were  first  conducted  to  justify  that,  in  the  case  of  a  homogeneous  half-space, 
the  grid  dispersion  only  modulates  the  spectra  of  the  propagating  Rayleigh  wave  (Figure  2). 
No  significant  attenuation  of  the  Rayleigh  wave  could  be  measured  for  the  half-space.  Follow¬ 
ing  this  check.,  we  introduced  the  simple  ramp.  We  measured  the  transmission  and  reflection 
of  Rayleigh  waves  incident  upon  simple  uphill  or  downhill  slopes  with  variable  height.  The 
reflection  and  transmission  coefficients  are  dependent  on  the  ratio  of  the  ramp  height  to  the 
wavelength  as  well  as  the  direction  of  the  ramp  as  shown  in  Table  1  for  constant  slope.  For  the 
same  height,  the  downgoing  slope  always  causes  more  reflection  and  less  transmission  as  com¬ 
pared  to  the  upgoing  slope.  Reflection  is  always  less  than  transmission.  For  s'opes  of  the 
same  direction,  larger  scale  variations  in  topography  produce  smaller  transmission  coefficients. 


Here  we  define  the  Rayleigh  wave  energy  transmission  coefficient,  T,  as  T(f)  = 


Power, raJf) 
Power 


The  reflection  coefficient,  R,  is  defined  in  an  analogous  way:  R(f)  = 


Powerreftf) 
r™’ermc(f)  ' 


Power  spec¬ 


tra  of  ihe  strain  (  dilatation  and  rotation)  were  used  to  compute  the  coefficients  so  as  not  to  be 
biased  by  the  dispersion  effects.  The  strain  in  each  case  was  "measured"  one  grid  point  below 
the  flat  surface  away  from  the  slopes  in  question. 
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TABLE  1.  Energy  Transmission  &  Reflection  Coefficients 

h/A 

uphill* 

downhill* 

T(lHz) 

R(lHz) 

1-T-R 

T(lHz) 

R(lHz) 

1-T-R 

1/16 

.87 

.01 

.12 

.75 

.016 

.22 

1/4 

.51 

.003 

.48 

.44 

.008 

.55 

1/2 

.54 

.004 

.45 

.16 

.009 

.83 

1 

.62 

.0036 

.38 

.06 

.035 

.93 

2 

.22 

— 

.0033 

.88 

.034 

.015 

.95 

*all  ramps  have  45°  slope. 


A  significant  conclusion  is  immediately  apparent  from  Table  1,  that  the  reflected  energy 
and  transmitted  energy  is  far  from  100%  of  the  incident  energy  and  that  bodywave  conversion 
must  be  very  efficient.  Clearly  55%  of  the  Rayleigh  wave  energy  is  lost  to  bodywave  conver¬ 
sion  by  passage  down  a  A/4  slope  of  45°. 

Figures  3  through  6  show  snapshots  of  the  displacement  field  due  to  a  Rayleigh  wave 
incident  upon  45°  ramps  of  various  heights.  Successive  frames  are  separated  by  intervals  of  2 
seconds  and  the  incident  Rayleigh  wave  always  travels  from  left  to  right.  Displacements  are 
proportional  to  the  darkness  of  the  plot  and  are  normalized  to  the  maximum  in  each  frame. 
The  diffraction  pattern  becomes  much  more  complicated  when  the  the  height  of  the  ramp  is 
larger.  This  is  obvious  in  comparing  Figure  3  (h  =  0.5A  )  and  Figure  5  (h  =  2A  ).  In  all  cases 
the  corners  act  as  point  sources  radiating  converted  bodywaves  (mainly  S  wave).  When  the 
incident  Rayleigh  wave  encounters  a  rising  hill,  some  energy  is  converted  to  P  wave  energy 
and  propagates  in  essentially  the  same  direction  as  the  incident  Rayleigh  wave  (Figure  6  and 
8),  while  in  the  case  of  going  down  hill,  the  converted  P  wave  scatters  in  all  possible  direc¬ 
tions  (Figures  3,  4,  5,  and  7).  In  the  case  of  the  2A  downhill  slope,  the  amplitude  spectra! 
ratios  of  the  forward  scattered  (with  30°  takeoff  angle)  P-  and  S-  waves  over  the  backward 
scattered  (with  -30°  takeoff  angle)  P-  and  S-  waves  near  1  Hz  are  0.73  and  2.82  respectively. 
This  is  to  say  that  the  forward  scattered  P  waves  are  0.73  times  as  large  as  the  back  scattered 
P  waves.  The  corresponding  measurements  in  the  case  of  going  uphill  are  1.4  and  1.3  for  P 
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and  S  waves  respectively.  All  these  observations  indicate  that  Rayleigh-to-bodywave  conver¬ 
sion  is  directional  dependent.  Since  Lg  is  dominated  by  S  waves  with  a  slowness  correspond¬ 
ing  to  a  takeoff  angle  in  our  simulation  of  about  30°,  we  see  that  the  Rayleigh-to-Lg  conver¬ 
sion  is  preferentially  in  the  forward  direction  by  an  average  factor  of  2  in  amplitude  (in  this  h 
=  2X  case). 

Figure  7  and  8  show  the  snapshots  of  a  Rayleigh  wave  incident  upon  a  valley  and  a  tra¬ 
pezoidal  mesa  of  depth  IX  respectively.  The  scattering  is  dominated  by  the  first  corner  or 
ramp  encountered  during  the  propagation  and  hence  the  diffraction  pattern  of  the  valley  is  very 
similar  to  that  of  single  downgoing  ramp  (Figure  4).  Likewise,  the  diffraction  pattern  of  the 
mesa  is  similar  to  the  upgoing  ramp  (Figure  6)  as  expected. 

The  transmitted  and  reflected  Rayleigh  waves  can  be  seen  in  Figures  9A,B  and  10A.B  for 
the  2X  uphill  and  downhill  cases.  The  amplitudes  of  these  waves  are  in  close  agreement  with 
the  reflection  and  transmission  coefficients  determined  from  the  power  spectra. 

When  we  examine  the  spectra  of  the  forward-  and  back-scattered  P  and  S  waves,  we  find 
that  we  can  separate  the  high-  and  low-frequency  scattering  regimes  as,  X  <  h,  and  X  >  h 
respectively.  From  these  simple  observations  we  can  compare  the  forward-  and  back-scattered 
P  and  S  wavefields  as  a  function  of  h/X.  The  converted  P  waves  are  predominately  forward 
scatter  (about  10  times  larger  than  the  back  scattered  P  waves)  for  X  <  h,  and  more  isotropic  at 
low  frequencies,  for  X  >  h.  This  applies  to  Rayleigh  waves  converging  on  an  uphill  or 
downhill  slope.  In  contrast,  the  converted  S  waves  are  roughly  the  same  amplitude  in  the  for¬ 
ward  and  backward  directions  for  both  high  and  low  frequencies  for  Rayleigh  waves  incident 
on  a  downhill  slope.  The  converted  S  waves  from  uphill  slopes  are  dominated  by  forward- 
scattered  waves  at  high  frequencies,  and  dominated  by  back -scattered  waves  at  low  frequen 
cies. 
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If  we  examine  the  P  and  S  waves  radiated  by  the  ramps  in  Figures  3  through  6  and  Fig¬ 
ures  11  through  14.  we  see  that  there  is  an  SV  node  in  the  vertical  direction  from  each 
diffractor  while  the  P  wave  radiation  is  of  constant  polarity.  Intuitively,  this  follows  from  the 
retrograde  elliptical  P-SV  motion  of  the  incident  Rayleigh  wave.  The  fundamental  mode  has 
no  shear  strain  component  on  a  horizontal  plane  in  the  halfspace.  This  fundamentally 
influences  the  radiation  patterns  of  the  Rayleigh-to-body  wave  diffractors.  The  diffraction  pat¬ 
terns  resemble  those  of  a  point  force  directed  in  a  direction  bisecting  the  local  normals  at  the 
corner.  More  SV  energy  will  therefore  tend  to  be  trapped  in  the  crustal  waveguide  and  the  P 
waves  are  more  likely  to  escape  the  crustal  waveguide  at  steeper  angles  of  incidence. 

Wavenumber  spectra  from  the  seismic  sections  of  Figures  13  and  14  are  shown  for  1,  2, 
3,  and  4  Hz  in  Figures  15A,B  and  16A.B.  They  indicate  that  most  of  the  energy  of  both  the  P 
and  S  waves  is  steeply  incident.  The  principle  P  and  S  waves  go  from  steeply  incident  in  the 
backward  direction  at  1  Hz  to  steeply  incident  in  th  forward  direction  as  frequency  increases 
from  1  to  4  Hz. 

These  experiments  with  simple  topographic  profiles  are  instructive  since  they  show  us 
that  Rayleigh-to-bodywave  conversion  occurs  near  corners  (changes  in  slope),  is  frequency 
dependent,  directionally  dependent,  and  in  general  depends  on  the  sign  in  the  change  of  slope. 
These  phenomena  are  more  difficult  to  observe  in  more  complicated  topographic  profiles. 
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NUMERICAL  EXPERIMENTS:  ROUGH  TOPOGRAPHIC  PROFILES 

In  order  to  generate  random  topography,  a  random  number  generator  is  used  to  produce  a 
Markov  sequence  of  values  0,  1,  or  -1  with  a  chosen  transition  probability  matrix.  The  associ¬ 
ated  random  walk  (integrated  sequence)  is  then  derived.  The  finite-difference  scheme  must 
discretize  the  topography  into  polygons,  therefore  the  smallest  scale  length  of  the  topography  is 
always  twice  the  horizontal  grid  spacing.  The  largest  scale  length  is  limited  by  the  size  of  the 
grid  used  in  the  simulation.  Several  typical  topographic  profiles  derived  by  this  method  are 
shown  in  Figures  17A  through  17E. 

15  rough  topographic  profiles  were  inserted  into  the  grid  with  flat  ground  concatenated  on 
both  sides  to  simulate  the  propagation  of  Rayleigh  waves  across  the  grid.  The  initial  incident 
waveform,  the  reflected  and  the  transmitted  waveforms  were  then  measured  at  near  surface 
"sensors"  away  from  the  rough  topography.  The  grids  are  26km  in  depth,  with  grid  spacing  of 
0.1  km  and  widths  between  34  to  39  km  depending  on  the  profile.  Transmission  coefficients 
and  the  attenuation  factor  Q(f)  were  then  computed  in  the  frequency  domain.  The  results  are 
shown  in  Figure  17A  through  17E.  The  useful  bandwidth  is  limited  to  frequencies  below  4  Hz 
due  to  the  limitation  of  2nd  order  LFD  scheme.  Although  we  quote  dimensions  in  km  and 
velocities  in  km/s,  it  is  possible  to  scale  these  simulations  to  higher  or  lower  frequencies  as 
desired. 


Table  2  gives  the  statistics  of  several  topographic  profiles  used  and  the  corresponding 
attenuation  factor  Q  at  the  dominant  frequency.  Here  we  computed  the  attenuation  factor  Q 

i  VK 

via  Q  (f)  =  -ln(T(f))( - ).  The  roughness  coefficient  Qj  (  or  normalized  total  variation)  is 


defined  by  12)  s - - -  where  //,-  is  the  elevation  at  /'th  column  of  the  nontrivial  por- 
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tion  of  the  topography,  X  the  traveling  distance.  This  statistic  reflects  the  idea  that  in  a  smooth 
sequence  of  points  the  differences  between  successive  values  are  all  small,  and  hence  the 
numerator  will  be  small.  The  denominator  is  simply  a  scaling  factor.  Further  investigation  on 
whether  such  £2t,  as  a  measure  of  roughness,  can  characterize  all  general  rough  topographic 
features  is  necessary. 


In  these  numerical  tests,  smoother  topographic  profiles  have  larger  Q  as  expected.  For 
the  specific  profiles  listed  in  Table  2,  the  ordering  is  preserved  if  Q!  were  replaced  by  the  nor- 

ft, _ // 

malized  l2  norm  of  the  elevation  change,  which  is  defined  by  - - — — — .  Q  esti¬ 


mates  are  much  the  same  whether  shear  or  dilatadonal  strain  energy  were  used  to  measure 
attenuation.  Several  calculations  were  run  at  grid  sizes  differing  by  factors  of  2  and  1.428. 
The  resultant  Q’s  were  nearly  identical. 


TABLE  2.  Rough  Topo 

graphic  Profiles 

topography 

span(km) 

Q(.8Hz)* 

3 

8.7 

.18 

629 

15 

11.4 

.39 

184 

5 

11.9 

.40 

37 

4 

12.1 

.46 

33 

10 

11.5 

.49 

28 

6 

12.2 

.52 

25 

11 

11.2 

.54 

21 

8 

12.0 

.58 

15 

*  ( VR  =2.71  km/s,  A,  =3. 2km) 


A  horizontal  array  of  64  sensors  are  equally  spaced  near  the  bottom  of  the  grid  to  record 
the  converted  P  (dilatation)  and  SV  (rotation)  wavefields.  Strain  energy  ratios  were  computed 
for  the  transmitted  Rayleigh  wave  to  the  incident  Rayleigh  wave  in  order  to  estimate  the 
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energy  loss.  Power  spectra  were  used  in  order  to  avoid  the  dispersive  effects  of  the  scattering. 
Power  ratio  of  the  converted  P  (S)  waves  to  the  incident  Rayleigh  wave  are  used  as  a  measure 
of  scattering.  Furthermore,  a  frequency-wavenumber  analysis  technique  is  used  to  produce  F- 
K  plots  of  the  scattered  dilatational  and  rotational  bodywavefield  as  observed  near  the  bottom 
of  the  grid. 

Seismic  sections  from  a  linear  array  ( Ax  -  1  km)  at  depth  within  the  grids  for  profiles  #5, 
#10,  and  #15  are  shown  in  Figures  18,  19  and  20.  The  P  and  S  wavefields  can  be  seen  to  be 
steeply  incident  upon  the  array  from  numerous  point  diffractors  in  the  topographic  profile.  The 
coda  lasts  as  long  as  it  takes  the  Rayleigh  wave  to  propagate  across  the  profile.  Wavenumber 
spectra  from  a  linear  array  (Ax  =  0.3  km)  of  sensors  near  the  bottom  of  the  grid  are  shown  for 
profiles  #5,  #10,  and  #15  in  Figures  21,  22,  and  23.  At  1  Hz,  the  wavenumber  spectra  show  a 
broad  pronounced  maximum  corresponding  to  nearly  horizontally  propagating  P  or  S  waves. 
The  1  Hz  P  and  S  waves  are  predominately  back-scattered  waves.  At  2  Hz,  the  maximum 
broadens  considerably  and  the  forward-scattered  waves  are  roughly  equal  in  power  to  the 
back-scattered  waves.  The  S  wave  spectra  tend  to  show  a  notch  near  vertical  incidence.  At  3 
and  4  Hz,  the  body  waves  are  mostly  forward-scattered  waves. 


The  spectra  of  the  bodywaves  reflect  the  spatial  attenuation  of  the  Rayleigh  wave.  The 
strain  energy  partitioning  between  the  P  and  S  waves  is  approximately  l-to-10.  when  the  areas 
under  the  wavenumber  power  spectra  are  accounted  for.  This  is  roughly  independent  of  fre¬ 
quency.  The  amplitude  ratio  between  the  S  and  P  waves  is  therefore  roughly  consistent  with  a 


ratio  of  (— ):  as  w°uld  be  expected  from  a  distribution  of  point  forces  rather  than  a  ratio  of 
(-j^-)3  as  would  be  expected  from  a  distribution  of  shear  stresses  at  the  surface.  Recall  that  the 


P-to-S  ratios  from  the  simple  45°  ramps  were  much  smaller  and  had  strong  frequency  depen¬ 


dence 


Scientific  Report  3 


12 


NOVEMBER  1986 


AFGL-TR-86-0269 


R-Wave  Scattering  by  Topography 


TGAL-86-9 


The  relationship  between  the  Q~l(f)  and  the  wavenumber  spectrum  of  the  topography  is 
not  as  clear  as  we  would  like.  Table  2  shows  a  subset  of  C21  versus  Q  for  f  =  0.8  Hz  derived 
from  several  topographic  profiles.  However,  another  subset  of  random  topographies  were 
found  to  have  a  different  relationship  between  the  statistical  moments  and  the  Q  at  0.8  Hz. 
We  have  plotted  in  Figure  10,  the  topographies,  //( x),  their  wavenumber  spectra,  H(k),  and  the 
estimated  Q(f)  for  15  different  topographic  realizations.  Some  of  the  statistics  of  these  topo¬ 
graphies  are  listed  in  Table  3  along  with  the  Q’s  at  1,  2  and  3  Hz  where  Q,  and  fi2  are 

defined  as  before,  and  Q3  = - — —  with  H  the  average  elevation  of  the  topography.  Q, 

I(//,-//)2 

is  also  called  Durbin-Watson  statistic  or  von  Neumann  ratio  (c.f  Bloomfield,  1975). 


TABLE  3.  Rough  Topographic  Profiles 


graphy 

n. 

Q(lHz) 

Q(2Hz) 

1 

.40 

.06 

.0063 

107 

74 

2 

.25 

.04 

.16 

75 

54 

3 

.18 

.05 

.15 

630 

63 

4 

.46 

.06 

.012 

44 

33 

5 

.40 

.06 

.029 

118 

94 

6 

.52 

.07 

.024 

23 

27 

7 

.47 

.06 

.0084 

19 

27 

8 

.58 

.07 

.012 

13 

31 

9 

.42 

.06 

.0096 

24 

36 

10 

.49 

.07 

.013 

30 

36 

11 

.54 

.07 

.012 

23 

36 

12 

.41 

.06 

.08 

211 

71 

13 

.41 

.06 

.08 

100 

57 

14 

. 

.42 

.06 

.023 

97 

80 

15 

.39 

.06 

.041 

104 

49 
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DISCUSSION  AND  CONCLUSIONS 

We  see  that  moderate  topography  on  a  half  space  can  substantially  attenuate  a  short- 
period  fundamental  Rayleigh  wave.  A  substantial  portion  of  the  energy  is  lost  to  bodywave 
conversion.  We  can  expect  to  see  P  and  S  waves  radiated  as  coda  to  teleseismic  distances  by 
short-period  fundamental  Rayleigh  waves  and  some  P-SV  energy  trapped  in  the  crustal 
waveguide  as  Lg.  A  significant  result  of  the  topographic  simulations  is  that  reflection  of  Ray¬ 
leigh  waves  by  topography  at  normal  incidence  is  very  inefficient.  Therefore  we  should  not 
expect  to  see  Rayleigh  back-scattering  as  a  significant  contributor  to  the  multiple  scattering  of 
fundamental  Rayleigh  waves  that  populate  P  coda. 

We  have  performed  the  2-D  numerical  Rayleigh-wave  scattering  experiments  both  upon 
simple  ramps  as  well  as  more  complicated  random  topographies.  The  simple  ramps  show  that 
reflection  is  an  inefficient  process  and  that  the  bulk  of  the  energy  that  is  not  transmitted  across 
a  step  is  converted  to  bodywaves.  The  transmission  through  a  step  is  frequency  dependent  and 
depends  on  the  direction  of  the  step.  The  bodywaves  appear  to  radiate  from  the  changes  in 
slope  and  to  first  order  can  be  visualized  as  from  a  point  force  diffractor.  The  S  waves  have  a 
null  in  the  direction  of  the  local  normal  similar  to  the  radiation  pattern  of  a  vertical  point  force 
At  low  frequencies  (h  <  X),  the  bodywaves  are  predominately  back-scattered,  while  at  high  fre¬ 
quencies  (h  >  X),  there  is  a  predominance  of  forward-scattering. 

More  complicated  topographic  profiles  contain  many  scale  lengths  and  introduce  broad¬ 
band  effects.  The  Rayleigh  wave  attenuation  is  a  complicated  function  frequency  but  is  a  max¬ 
imum  in  the  range  where  the  wavenumber  spectra  of  the  topography  matches  the  wavelength 
of  the  incident  Rayleigh  wave. 
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The  uncertainties  in  the  Q  estimates  are  often  approaching  100%  for  the  Q~x  values  less 
than  0.01  in  Figure  17.  A  expected,  the  lowest  Q’s  are  always  the  best  determined.  When 
several  topographic  profiles  are  selected  from  an  ensemble  with  the  same  statistics  we  see  that 
the  Q  also  has  a  characteristic  variance.  For  example,  the  four  profiles  #12,  #13,  #14,  and  #15 
were  generated  by  using  different  seeds  and  have  similar  spectra  and  Q(f)  with  an  average  Q 
of  90  in  the  1  to  2  Hz  band  with  a  standard  deviation  of  30%. 

The  dilatational  and  rotational  wavenumber  spectra  indicate  that  the  scattered  P  and  S 
waves  are  scattered  in  the  backward  direction  at  low  frequencies  (1  Hz)  and  in  the  forward 
direction  at  high  frequencies  (4  Hz).  Since  the  profiles  are  band-limited,  X  <  total  span  and  X 
>  0.2  km,  the  intermediate  frequencies  (1  to  3  Hz)  are  probably  the  most  reliable  (0.9  <  X  < 
2.7).  In  this  bandwidth,  we  see  that  the  ratio  of  the  P-to-S  waves  is  nearly  constant  with  a 


Ot  o 

value  near  (— )  as  if  the  apparent  source  were  a  distribution  of  point  forces,  and  the  scattering 

r 

is  a  mixture  of  forward-  and  back-scattering.  While  apparently  more  SV  energy  may  be 
trapped  in  the  crustal  waveguide,  the  clearly  defined  node  in  the  vertical  direction  is  less  well 


defined  for  random  rough  topography  and  the  SV  energy  is  more  evenly  distributed  in  takeoff 


angle  or  slowness. 


The  results  presented  offer  a  beginning  approach  to  effects  of  scattering  by  topography 
and  heterogeneity  upon  surface  waves.  We  have  assumed  that  topography  can  be  approximated 
by  level  elevations  separated  by  sections  of  45°  slope.  Our  principle  limitation  is  that  we 
require  segments  of  topography  to  exceed  2Ax  (Jih  et  al ,  1986)  and  that  Ax  must  be  less  than 
XJ 30  (Fuyuki  and  Matsumoto.1980).  A  complete  exploration  of  the  problem  will  require  varia¬ 
tion  of  the  P  and  S  velocities,  near  surface  velocity  gradients,  as  well  as  crustal  velocity 
heterogeneity.  The  short-period  Rayleigh  will  not  efficiently  propagate  across  rough  topography 
and  reflection  is  small  Unless  other  mechanisms  are  found  to  tr  p  P-SV  energy  within  the 
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crust,  the  bulk  of  the  energy  will  leak  out  of  the  waveguide.  This  can  be  mitigated  by  lowering 
the  near  surface  velocity  to  spread  the  slowness  space  available  to  the  converted  P  and  SV 

waves,  increasing  the  ratio,  introduction  of  gradients  near  the  surface  to  create  higher  order 

modes,  and/or  velocity  heterogeneity  to  scatter  P  and  SV  energy  back  into  the  waveguide.  In 
short,  all  these  mechanisms  can  only  increase  the  Rayleigh-to-Lg  coupling. 


i 
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FIGURE  CAPTIONS 

Figure  1A  Wavenumber  spectrum  S(A :)  =  {idk^expl  I — (^/Aq)2)  of  Ricker  wavelet  where  k  is 
wavenumber  and  Icq  =  1 .963  radians/km  the  maximum  amplitude  wavenumber. 

Figure  IB  Ricker  wavelet  in  the  spatial  domain.  IFFT  was  applied  on  S(k)  to  speed  up  the 
wavelet  generation. 

Figure  1C  Rayleigh  wave  displacement  eigenfunctions  for  the  vertical  and  horizontal  com¬ 
ponents  of  motion. 

Figure  2  Spectral  ratio  of  the  transmitted  Rayleigh  waves  over  the  incident  Rayleigh  waves  in 
a  homogeneous  half-space.  Grid  dispersion  only  modulates  the  spectra. 

Figure  3  The  snapshots  of  the  displacement  field  due  to  Rayleigh  wave  incident  on  a  45- 
degree  ramp  of  height  h  =  0.5  X.  Successive  frames,  A  through  F,  are  separated  by  2  sec 
intervals.  Displacements  are  proportional  to  the  darkness  of  the  plot  and  are  normalized  to  the 
maximum  in  each  frame.  Note  the  P  and  S  waves  generated  by  the  two  comers  of  the  ramp. 

Figure  4  Same  as  Figure  3  except  h  =  X.  The  diffracted  wavefield  is  more  complicated  as 
compared  to  the  case  h  =  0.5  X  (Figure  3).  There  is  a  minimum  in  the  scattered  S  wave  in 
the  near  vertical  direction. 

Figure  5  Same  as  Figure  3  except  h  =  2  X.  The  diffracted  wavefield  is  more  complicated  as 
compared  to  the  cases  of  smaller  h  (Figure  3  and  4)  or  the  case  of  upgoing  cliff  with  the  same 
height  (Figure  6).  Note  the  constructive  interference  of  the  forward-scattered  S  wave  by  the 
two  corners  and  destructive  interference  for  the  back- scattered  S  wave. 

Figure  6  Same  as  Figure  5  except  the  direction  of  the  ramp  is  reversed.  Note  that  the 
diffracted  P  and  S  waves  propagate  mainly  to  the  right  while  for  downgoing  ramps  (Figure  3, 
4,  5)  or  valley  (Figure  7)  the  diffracted  bodywaves  propagate  more  isotropically. 

Figure  7  Rayleigh  wave  incident  upon  a  valley  of  depth  1  X.  The  scattering  is  dominated  by 
the  first  downgoing  ramp  encountered  during  the  propagation  and  hence  the  diffraction  pattern 
is  very  similar  to  Figure  4. 

Figure  8  Rayleigh  wave  incident  upon  a  trapezoidal  mesa  of  height  1  X.  The  high-frequency 
fraction  of  the  Rayleigh  wave  propagates  along  the  mesa  surface  while  the  low-frequency  frac¬ 
tion  tunnels  through  the  mesa. 

Figure  9A  Synthetic  near-surface  vertical  displacements  (upper)  and  horizontal  displacements 
(below)  for  a  Rayleigh  wave  incident  on  a  2X  deep  downhill.  Location  of  sensors  are  shown 
in  Figure  5A.  Note  the  phase  dispersion  of  the  transmitted  Rayleigh  wave.  Two  reflected 
waves  are  seen  from  each  of  the  two  corners. 

Figure  9B  Synthetic  near-surface  dilatational  strain  (upper)  and  rotational  strain  (below)  for  a 
Rayleigh  wave  incident  on  a  2X  deep  downhill.  Location  of  sensors  are  shown  in  Figure  5A. 

Figure  10A  Same  as  9A  except  for  opposite  direction  of  the  ramp.  Location  of  sensors  are 
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shown  in  Figure  6A. 

Figure  10B  Same  as  9B  except  for  opposite  direction  of  the  ramp,  location  of  sensors  are 
shown  in  Figure  6A. 

Figure  11  Seismic  sections  recording  the  converted  P  wave  (dilatational  strain,  upper)  and  S 
wave  (rotational  strain,  lower)  at  a  line  of  32  sensors  near  the  bottom  of  the  grid  spaced  1  km 
apart  for  the  case  of  0.5X  deep  downhill.  See  Figure  3. 

Figure  12  Same  as  1 1  except  for  opposite  direction  of  the  ramp.  See  Figure  4. 

Figure  13  Same  as  1 1  except  for  2X  deep  downhill.  See  Figure  5. 

Figure  14  Same  as  12  except  for  opposite  direction  of  the  ramp.  See  Figure  6. 

Figure  15A  Log  power  horizontal  wavenumber  spectra  for  the  dilatational  wavefield  recorded 
across  the  bortom  linear  array  for  the  2X  downgoing  ramp  case.  Power  spectra  at  1,2,  3,  4Hz 
are  normalized  to  the  maximum.  Normally  incident  waves  have  a  horizontal  wavenumber  of 
zero  Wavenumbers  are  indicated  for  horizontally  propagating  P  waves  Forward -scattered 
waves  have  negative  wavenumbers,  back -scattered  waves  have  positive  wavenumbers.  Note 
that  the  dilatational  power  is  confined  to  steep  angles  of  incidence  across  the  array,  slightly 
back-scattered  at  1  Hz.  and  more  isotropic  at  4  Hz. 

Figure  15B  la>g  power  wavenumber  spectra  for  the  rotational  wavefield  recorded  across  the 
bottom  of  the  grid  for  the  2X  downgoing  ramp  case.  S  waves  are  back -scattered  with  steep 
angle  of  incidence  at  1  H’  ind  more  isotropic  at  4  Hz  but  still  steeply  incident. 

f  igure  16A  Same  as  ICA  except  for  opposite  direction  of  the  ramp. 

Figure  16R  Same  as  l^B  except  for  opposite  direction  of  the  ramp. 

Figures  17A,B,(  ,I),F.  t  opography,  H(x),  on  the  left,  log-amplitude  wavenumber  spectra  of  the 
topography  in  the  center,  and  the  estimated  Rayleigh-wave  log(l/Q(f))  on  the  right  for  topo¬ 
graphic  profiles  #1  through  #15 

Figure  18  Seismic  sections  recording  the  converted  P  wave  (dilatational  strain,  upper)  and  S 
wave  (rotational  strain,  below)  at  a  line  of  32  sensors  located  near  the  bottom  of  the  grid 
spaced  1  km  apart  for  rough  topography  #5  Complicated  interference  patterns  are  evident  for 
both  P  and  S  wavetields 

Figure  19  Same  as  18  except  for  rough  topographic  profile  #10. 

Figure  20  Same  as  18  except  for  rough  topographic  profile  #15 

Figure  21 A  Frequency  wavenumber  spectra  for  the  dilatational  strain  field  (P-wave)  recorded 
near  the  bottom  of  the  grid  with  topography  #5  on  the  top  of  the  grid.  The  dilatational  strain 
energy  is  largely  confined  to  P  wave  slowness  across  the  array.  P  wave  energy  shifts  from 
back-scattered  to  forward-scattered  from  I  to  4  Hz. 
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Figure  21 B  Frequency  wavenumber  spectra  for  the  rotational  strain  field  (S-wave)  recorded 
near  the  bottom  of  the  grid  with  topography  #S  on  the  top.  The  rotational  strain  energy  is 
largely  confined  to  S  wave  slowness  across  the  array.  S  wave  energy  shifts  from  back- 
scattered  to  forward-scattered  from  1  to  4  Hz. 

Figure  22A  Same  as  21 A  except  for  rough  topographic  profile  #10. 

Figure  22B  Same  as  21 B  except  for  rough  topographic  profile  #10. 

Figure  23A  Same  as  21A  except  for  rough  topographic  profile  #15. 

Figure  23B  Same  as  21B  except  for  rough  topographic  profile  #15. 
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Figure  1C 
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NORMALIZED  Loq(POWER)  WAVENUMBER  SPECTRA 

Figure  15A 
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NORMAL  1  Z  L  D  LoqCPOWtR]  WAVENUMBER  SPECTRA 

Figure  16B 
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Figure  17B 
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NORMALIZED  Log  C POWER ;  WAVENUMBER  SPECIE 

Figure  21 A 
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